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Abstract 



A generic problem of high frequency wave propagation along a metallic strip 
in parallel above a PEC ground plane is considered. The wave is excited by an 
elemental electric dipole at an arbitrary location above the PEC plane. The full 
wave problem, for arbitrary widths of the strip, is solved by means of a mode 
matching approach and expansion of the strip surface current into Chebyshev poly- 
nomials. For narrow strips, an approximate method using only longitudinal currents 
is derived, and compared numerically with the full wave method. Utilizing the con- 
cept of equivalent radius, the approximate method for narrow strips is evaluated 
numerically against results for thin circular wires. It is concluded that the approxi- 
mate method is suitable for handling multiple wires in layered structures, wherefore 
the method has potential usefulness for estimating long range propagation of high 
frequency waves in wire structures like power lines and railway feeding systems, 
containing over- head wires and wires submerged into ground. 

1 Introduction 

The analysis of electromagnetic waves propagating along wires parallell with interfaces, 
separating regions having different electromagnetic properties, is of importance when 
considering electrical systems, like power lines submerged cables [HITT] and the power 
supply for railways [1]. In such power supply systems, the waves can be unintentional, 
triggered by internal and/or external disturbances [HISIITO], or intentional, for e.g. moni- 
toring and communication purposes [2llll|5]. 

The literature on methods for analyzing parallell conductor structures is vast, wherefore 
we will not give a comprehensive survey here; a selection of references can be found in 
the monograph [25] . In a low frequency application, it is suitable to model the structure 
as a multiconductor transmission line (MTL), and solve the propagation problem by 
means of the quasi-TEM mode theory [22l[23], generalized to account for the coupling to 
external sources [25] (chapter 1). The basic assumption is that at the highest frequency, 
the wavelength is large compared with the transverse dimensions of the system. For power 
lines and railway feeding systems one can thus reliably use MTL theory at frequencies 
below about 1-10 MHz. For frequencies at which the wavelength becomes comparable with 
the wire separations the presence of higher order radiating modes cannot be neglected. 
In such cases one can use full wave, or nearly full wave, methods, like generalized MTL 
methods [25] (chapter 4) or other methods, like the Finite Difference Time-Domain method 
[9] and integral equation based methods [6l[7]. 

Typically, the full wave methods are based on numerical discretization of wires of finite 
length, whereby for very long wires the required numerical accuracy may yield equations 
that are too massive from the computational point of view. In particular, this will be 
the case if one wants to study electromagnetic interference (EMI) and communication at 
frequencies of orders from 100 MHz to several GHz, where a typical power line becomes 
electrically very large also in the transversal directions. Hence, in order to estimate 
e.g. possible distances for GHz communication along wire structures and EMI effects 
at the same frequencies one needs methods capable of handling very long (in terms of 
wavelengths) wire structures. 
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One example of a nearly full wave method for an infinitely long uniform MTL is presented 
in [13] . The MTL is located above a lossy halfspace and illuminated by a plane wave with 
an arbitrary direction of incidence; remotely located dipole sources are also considered. 
The uniformity facilitates the use of a spatial Fourier transform whereby the problem 
can be solved in the transform space, i.e. the spectral domain. Except for the thin wire 
approximation, this method is exact. 

Spectral domain methods have been used frequently for studying propagation along MTL 
structures in the form of thin planar strip conductors on substrates, see e.g. Sect. 4.8 
in [26] for a calculation of the transmission line parameters for quasi-TEM modes propa- 
gating along a microstrip line. Hence, it is of interest to investigate whether these methods 
can be useful also for the above mentioned power line problems. A motivation for con- 
sidering methods developed for planar conductors is that if the thin wire approximation 
holds, and proximity effects can be neglected, a wire of quite arbitrary cross-section can 
be replaced with an equivalent strip conductor. For example, for a circular wire the 
equivalent strip has a width that is twice the diameter of the wire [21]. In [15], the two- 
dimensional scattering problem for a strip illuminated by a plane wave has been solved, for 
both polarizations of the incident electric field. Excitation/scattering problems for strip 
gratings using localized nearby sources have been considered in [19l[20] . In the context of 
power/railway lines, localized sources are of relevance when studying EMI emanating from 
e.g. spark-overs in insulators or from pantograph arcing [27]. Another situation where 
localized sources are of relevance is high frequency communication, where the guided wave 
on the wire structure can be excited by coupling to a suitably placed antenna. 

In this paper, we consider a metallic strip placed in parallel above a perfectly electrically 
conducting (PEC) plane. The configuration is excited by an electric dipole at an arbitrary 
location above the PEC plane. For the analysis, we use a spectral domain method, based 
on Fourier transforms in both of the horisontal directions. The analysis is systemized 
in terms of classical waveguide theory, where the spectral components are treated as 
waveguide modes in the vertical direction [2B] (Section 3.2b). Hence, the strip can be 
considered as a generalized diaphragm in a plane connecting two waveguide sections; a 
similar point of view has also been used in [19]. Like in [15], we will in the context of the 
strip model make no approximations, whereby the analysis holds for arbitrary widths of 
the strip, in terms of wavelengths. Furthermore, there are no restrictions on the separation 
distance between the wire and the PEC plane. 

The paper is organized as follows. In Section [2l we formulate the problem and derive the 
spectral method. For clarity, the derivation is carried out in some detail, with additional 
details in the appendices. The result is a linear system for determining the coefficients 
of the spectral surface current density on the strip. In Section [3], we present numerical 
results. First, we consider a wide strip, whereby we have to use the exact (truly full 
wave) method. Then we consider a narrow strip and investigate the applicability of 
certain approximations for that case. Finally, we evaluate the equivalent radius concept, 
by comparing the results for narrow strips with results for a circular wire: analytical 
results for the TEM-mode approximation and full-wave results obtained by the MoM 
based Numerical Electromagnetics Code (NEC) [31]. The method and the results are 
summarized and discussed in Section HJ 
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2 Problem formulation and theory 



The problem geometry is depicted in Figure [H There is a PEC surface in the plane z = 0. 
In the plane z = a there is an infinitely thin PEC strip that is infinite in the x-direction 
and with edges aX y = ±h. At the location ro = xqX + y^z + zqZ there is an electric 
dipole source with dipole moment p. The medium in the region 2 > is homogeneous, 
isotropic and in general lossy, described by the complex permittivity e and the complex 
permeability /i, which are in general frequency dependent. 

The problem is to find the fields and the strip surface current that results from to the 
radiation from the dipole. 
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2.1 Preliminaries 

We start from the Maxwell equations in a homogeneous isotropic medium: 

V X E{r) = -ikr]H{r) , (1) 
r/V X H{r) = +ikE{r) + r]J{r) , (2) 



where E and H are the electric and magnetic fields, respectively, k = Uy/ejl is the 
medium wave number, rj = \fjije is the medium wave impedance and J is the source 
current density, uj is the angular frequency, in a suppressed time-dependence. 

We decompose the radius vector into r = p + zz, where the transversal part p = xx-\-yy, 
and introduce the following pair of Fourier-transforms: 



E{kt,z) = j E{r)e^'''-Pds, (3) 

E{r) = 2^ I Eik,,z)e-^''^-fd\, (4) 

defined analogously for H, J etc.. In (|3]), ds = dxdy and S denotes the entire xy-plane. 
The spectral (Fourier) variable kt = k^x + kyij is a transverse wave- vector. In (jlj), 
d^kt = dkxdky and JC denotes the entire k^ky-plane. In the following, we refer to E{r) 
as the spatial electric field and to E{kt,z) as the spectral electric field, and similarly for 
other quantities. 
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2.2 TM- & TE-modes 



In any slab region, Zi < z < Z2, where J{r) = 0, it follows from ([T]) and ([2]) that the 
fields satisfy the homogeneous Helmholtz equation 



-E{r) 




0" 


H{r)_ 








(5) 



Using (j4]) to express E{r) and H{r) in ([5]), and using the property Ve j'^' '' = — jfet^ j^' '', 
it follows that the spectral fields satisfy 



-A; 



dz^ 





'E{k,,z) 
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H{k,,z)_ 








with solutions of the form 

E{k„z),H{k,,z)oce'^^'^' 
where the longitudinal wavenumber 



(6) 



(7) 



h = Vk^ - kl (8) 

Since any physical medium exhibits losses, it follows from the passivity conditions [30] that 
Im {e} < and that Im {/i} < 0, resulting in that k and kz both are complex. The branch 
of the square root in ([S]) is always chosen such that Im{/c^} < 0. Hence, e"-''^^^ represents 
a solution decaying in the +2;-direction, while e^^^'^ represents a solution decaying in the 
— z-direction. 

Using (j4]) to express the fields in Maxwell's equations ([1]) and ([2]), with J = and z- 
dependencies of the forms e^j'^^^, it is shown in Appendix |A] that the transversal spectral 
field components can be expressed in terms of the longitudinal components: 



[±k,KEf -kzx k,r]H^] 



= 7Y [±kzk,r]Hf + kz x k,Ef] 
kt 



(9) 
(10) 



where ± refers to the ±z-decaying solutions. 

From the structures of ^ and f[TU|) . we find it convenient to introduce the dimensionless 
functions 



9{kt) 



kz I 

kt 

fikt)Tz, 
^ - u 

~PyZ X /Ct, 

kt 



(11) 
(12) 

(13) 



subject to 



f^{k,)-g{k,) = f{k,)-g{k,) = 0. 



(14) 
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Hence, the transverse magnetic (TM) parts of the spectral fields are written 

r]H™ik,,z) = g{k,) [A+(fct)e-j'=^^ - ^-(fct) e+j'^^^] 
while the transverse electric (TE) parts of the spectral fields are written 
E'^^{k„z)=g{k,) [B+{k,)e-^'^' + B-{k,)e+^'''] 



(15) 
(16) 



(17) 
(18) 



vH^'^ih, z) = -B+ik,) /+(fct) e-j'^^^ + B-{k,) r{kt) e+^'^' 

For each value of k^ there is one TM-mode and one TE-mod^. A^{kt) and B^{kx) 
are mode coefficients, which are determined by the sources and the boundary conditions. 
Once the mode coefficients have been determined, the spatial fields are determined using 
(jlj), resulting in 



E{r) 
H(r) 



^ ^ [H™{k„ z) + /f™(fc„ ;.)] e-j^'-^d^A:,. 



(19) 
(20) 



2.3 The spectral fields from the dipole, the ground plane and 
the strip 

In Appendix [B] we show that the spectral fields from the dipole become 

E^''\k,, z) = A,{ku z) [H(z - zo) t{k,) + H(^o - z) r{k,)\ + fii(fct, z) g{h) (21) 

'i]H^'^\kt, z) = sgn{z - Zq) 

■ [A^{k„z)g{k,)-B^{k„z)[}l{z-Zo)f^ik,) + }^{zo-z)r{k,)]] (22) 

where H() denotes the Heaviside step-function, sgn() denotes the signum-function, and 
where for the principal directions of the dipole 



p = px 



p = py 



p = pz 




2e 
2ek, 

_iPh^kfPoQ-3k,\z-zo\ 

2e 

_]pkK^k,-po^-ik,_\z-zo\ 
2ek. 



(23) 



(24) 



iP^t Jkfp„ -}k^\z-ZQ\ 



sgn[z - zo) "e 



(25) 



^Whcn fct = there are no longitudinal fields, and instead two degenerate TEM-modcs (with = k). 



When enforcing the boundary conditions in a multilayered structure, it suffices to consider 
the transversal fields. For the dipole, the transversal fields become 

e[''\K z) = A,{ku z) f{kt) + fii(fct, z) 9{kt) (26) 

vHi'\k„ z) = sgn(^ - ^o) [Ai{k,, z) gik,) - Bi{k„ z) f{kt)] (27) 

The transversal fields originating from the charges and currents on the ground plane at 
z = decay away from the ground plane and are for z > thus written 

E['\h,z) = [Ao(fct)/(fct) + i?o(/2t)»(fct)]e-j'=^^ (28) 

TlHf\k,,z) = [Ao{h)gih) - i?o(fct) /(fct)]e^j'^^^ (29) 

The spectral transversal fields originating from the charges and currents on the strip in 
the plane z = a are written 

E[''\k„z) = [yl2(fct)/(fct) + fi2(fct)9(fct)]e-j'~~''""' (30) 

vH^ih, z) = sgn(z - a) [A^{k,) g{k,) - B^ih) /(fej] e^'^^l^-'^l (31) 

Note that since f{kt) ■ g{kt) = 0, f l30|) is the unique Ansatz to enforce the continuity of 
E^^\kt, z) across the plane z = a. 

2.4 Boundary conditions 

The free surface current density K{p) on the strip has the corresponding spectral surface 
current density 



(32) 



Applying the boundary condition H[^\kt, z = a~^) — H^^\kt, z = a ) = K{kt) x z, we 
obtain from (13T]) that 

2A2g - 2^2/ = i]K X z (33) 

Utihzing ([II]), ([I3D and ([HD, ([3SD yields 

_ g ■ {7]K X z) _ {zxg)- {r]K) _ k, ■ {r]K) _ k^K^ + kyKy 

V " 2^^ " 2^-'"^ 2k ^^^^ 

_ f-{r^Kxz) _ jzxf)- {r,K) _ {z x fcQ ■ (r/K) _ k^Ky - kyK, 

2f - 2P - 2k, -'"^ 2k. ^^^^ 

The condition Et = e'^'> + e'^^ + = at the plane z = 0, together with ffTI]) . fl2B]) . 
and f l5U]l . yields 

Ao = -Ai(z = 0)-A2e-j'^'^ (36) 

Bo = -B^iz = 0)-B2e-^^^'' (37) 
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Thus, using ([20]), (IMD-dSZD, (HH) and ([13]), we obtain that in the plane z = a 

the transversal spectral electric field becomes 



EUkt) = [Ai(fct,a) - Ai(fct,0)e-j'=^'^ + A2(fcO (l - e^^^^'^)] f{k 
+ [B,ih,a) - i?i(fct,0)e-j'=^'^ + B,ih) (l - e-j^^^'^)] g{k, 



^ _ g-j2feza-j kxkyKy{k t) (/c^ kl) Kx{kt 



k^k, (Ai(fct, g) - Ai(fct, 0) e-j^-") + kky {-Bijk,, a) + B,{k,, 0) e-^^-") >| , 

^ kl + kl r 

r _ g-j2fe,aj k^kyK^^kt) - (A:^ - fcg) i^^(fct) 

(Ai(fct, g) - Aijkt, 0) e-j^^-") + kk, (Ei(fct, g) - gi(fct, 0) e-^^-") „ 
^ kl + kl 

(3^ 



2.5 Spectral surface current densities 



In the plane z = a, the spatial surface current density K{x, y) is confined to the strip- 
conductor, with edges at ?/ = ±h. Explicitly taking into account the singular behaviour at 
the edges (see e.g. formula (4.13) in [33]), the components of the surface current density 
are expanded as 



Ky{x,y) 



oo 

H(/i-|y|)v/i-yV/^2^rfn(x) U„( 

m-\y\) 



(39) 



2^1 - yyh^ 



n=0 



n=0 



h 



y 
h 



(40) 



where T„() and U.„() denote the n:th order Chebyshev polynomials of first and second 
kind, respectively. {c„(x) ,(i„(x)}^Q are coefficient functions to be determined. 

Inserting (139]) and (140]) into (132!) and using formula 7.355 in [24], the components of the 
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spectral surface current density become 



Vr/i^ (-1)" [C2n{kx) hn{kyh) + j C2n+l{kx) hn+i{kyh)] 



(41) 



n=0 



oo 



^ J-h ^^l-yyh'' 



2 

7rh 



n=0 

oo 



= ^ ^ (-1)" d2n{kx) {hn{kyh) + hn+2{kyh)} 



n=0 



+ id2n+i{kx) {hn+i{kyh) + hn+sikyh)} 



oo 



n=0 



d2n{kx) 



(2n + 1) (fcj^/i) (2r2 + 2) J2„+2(fcy/i) 

— +]d2n+l[kx) — 

rvyit rvylv 



where the spectral coefficients 



/oo 
rf„(x)ej'^-^dx 
-oo 

and J„() denotes the n:th order Bessel function. 

2.6 Equation system for the spectral coefficients 

In the plane z = a, the spatial transversal electric field becomes 



(42) 



(43) 



Eta{x,y) 



47r2 



e-j'^^Mfc; / e-''^yydkyEta{k',,ky 



(44) 



On the strip, (— oo < x < oo) fl {—h < y < h)n{z = a), we have Eta{x, y) = 0. According 
to equation (4.9) in [33], it holds in the vicinity of an edge that oc y/r,Ey oc l/y/r, 
where r is the distance to the edge. Taking into account these properties, we enforce the 
vanishing tangential field on the strip through the following testing procedures: 

°° T (v/h) 

gjfc.x^^ / , "^^^ ' _ E^{x,y)dy = 0, -oo < k^ < oo, m = 0, 1, . . . 



oo 



H Vi - yVh' 



/oo ph 
eJ'^^Mx / ^/l - yyh^VUv/h) Ey{x, y) dy = 0, -oo < fc, < oo, 
-oo J —h 



m = 0, 1, 



(45) 



(46) 
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When inserting f l44p into f l45p and (146 p . the spatial integrals result in 

'h V 1 - yyh^ 

y^T~^VUy/h) e-^'yydy oc (m + 1) ^-I^±liMl 

-h ^yi^ 
and (H^ imply the following results in the spectral domain: 

/oo 
Ex{K, ky) Jmikyh) dky = 0, m = 0, 1, . . . 
'OO 



Hence, 



-oo 
f oo 



/oo 
Ey{kxi ky) 
-OO 



■Jm+1 

(kyh) 



dky = 0, m = 0, 1, . . . 



which, using (l38l) and (l23l)-(!25l). become 

/"°° _ g-j2fc^a-j (fe^ — kl) Kxjkt) — kxkyKy{ki) ^ 
J-oo 2kkz 



m {kyh) dky 



A;2 - k 



^ej'^-^^" / e^''yy°J^{kyh)dky< 



k, 

krpk' 



kz 



^ kx [sgn(a - ^o) e-j*^-l'^-^ol + e-j*^-('^+^")] 



(47) 
(48) 
(49) 



(50) 
(51) 



p = px 

p = py 

: p = pz 
(52) 



(1 - e-j2^^") 



(P - kl) Ky{k,) - kxkyKxih) 3„,+ ,{kyh) 



Iviih 



2er] 



I e^''yy'^J^+l{kyh)dky 



kyh 



-ikz\a~zo\ _ -jfcz(a+^o) 



] ■ p = px 



k^-k 

ky k^ h 



y ^Q-ikz\a-zo\ _ ^-ik^a+zo^ ■ p = py 



- [sgn \a - zo\ e-j'^^l'^-^ol + e-j'=^-('^+^'') 



] : p = pz 
(53) 



Next, we insert (HT]) and (H2|) into (152!) and (l53i) . and utilize the even and odd parities of 
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integer order Bessel functions, to obtain the following systems for {c2n{kx) , 'i2n+i(^x)}^o- 



-oo ^ _ g-j2fc^a 

7r/i2^(-l)" / ^ \ {k^ - kl) J2m{kyh) 32n{kyh) C2n{kx) 



n=Q 



kx 

j-^ (2n + 2) J27n{kyh) hn+2{kyh) d2n+i{kx) }dky 



er] 



hm{kyh) dky < 



X ^g-jfc,|a-.o| _ g-jfc.(a+.o)] cos(fey?/o) : P = pX 



] j sin(A;j^?/o) ■ P = py 



kx [sgn(a — Zq) e ■'''^1'' + e j'=^("+^o)j cos{kyyo) : p = pz 

(54) 



n=0 



Y fc2/^2 + 2) J2m.+2(/i;s//i) J2n+2(^y/i) C?2n+l(^:i 



■hm+2{kyh) J2n{kyh) C2n{kx) \dky 



er] 



kx 
kyh 



-ik^\a-zo\ _ -ik^(a+ZQ) 



] cos{kyyo) ■.p = px 



hm+2{kyh) dky < 



k^ - kl 



] j sin(A;y?/o) ■ P = py 



^ [sgn(a - zo) e-j^'^'"-""l + e-j'=^("+"»)] cos(A;y2/o) 



p = pz 

(55) 



and into the following systems for {c2n+i{kx) , (^2n(^2:)}^o- 



ixh\{-iy- / \]{k'^ - k"^^ :i2m+i{kyh):i2n+i{kyh)c2n+i{kx) 

n=0 "^0 

A; 

^ (2ra + 1) hm+iikyh) 

J2n+1 [kyh) d2n{kx) \dky 



h 



erj 



e - k 



'x ^g-jfc.|a-.'o| _ g 



-]kz{a+zo) 



] j sin(A;j^2/o) : p = px 



hm.+iikyh) dky < 



] cos{kyyo) ■.p = py 



kx [sgn(a - zo) + Q-]k.{a+zo)■^^ j sin(A;yyo) 



: p = pz 

(56) 
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oo 

nhY,i-l) 

n=0 



1 — Q-i'^k^a ( k 



^ {2n + 1) :i2m+i{kyh) hn+iikyh) 4n(fcx) 



er] 



i-^hm+l{kyh) J2n+l{kyh) C2n+l{kx) }dky 



-jfc^la-zol 



i^x r 

[e-^ ^ 



J2m+i{kyh) dky < 



jfc.(a+^o)]jsin(fc,yo) ■.P = PX 
y_ ^g-jfc.|a-2o| _ g-jfc4a+^o)] cos(/t3,?/o) : P = PS/ 



[sgn(a - zo) e-j*^^l"-""l + e-j'=^('^+^«)] j ^m.{kyy^) : p = pz 

(57) 



Summarizing, ( l5^ -( l571) are infinitely large linear equation systems for the spectral coeffi- 
cients Cn{kx) , dn{kx)- The elements in the system matrices and driving terms are obtained 
through various integrals with respect to the spectral variable ky, with the other spectral 
variable kr ClS db parameter. 

For the numerical solution, (!54|) - (!57|) are truncated into a finite number of coefficients 
and testing functions, and solved for different values of kx- With the spectral coefficients 
determined, the spectral surface current density follows from (14T]) and (H2|) . the mode 
coefficients from (15^ - (l?n) . and the spectral fields from (fT^ - (IT5l) . The details of the 
numerical evaluation of the spectral integrals are described in Appendix O Finally, the 
spatial fields follow from ffT9|) and f|20|). 

The coefficients for the spatial surface current density follow from the inverse Fourier 
transform of (l43l) : 



Cn{x) = — / Cn{kx) G ^ ^ dkx, dn{x) = — / dn{kx) G ^ ^ dkx 



(5^ 



With To = 1, it follows from ( l39l) and the orthogonality relation for the Chebyshev 
polynomials that the total current flowing along the strip becomes 



Kx) 



Kxix,y) dy 



dy 



oo 



T / \ To(m)T„,(m) ^ 

ri=0 "^"1 V 1 — M „_n 



(59) 



n=0 



3 Numerical results 



In this section, we evaluate the method through several numerical examples. In all ex- 
amples, the dipole source has the x-coordinate Xq = m, the frequency / = 300 MHz, 
and the medium parameters e = eg (1 ~ jlO~^) , = /iq (1 — jlO~^). Hence, it follows that 
the wavelength A ~ 1 m and that attenuation due to medium losses have minor impact 
within distances of the order of some hundred meter^. 

^In more realistic cases, attenuation is due to a finite ground conductivity that will mask the losses 
in the air. The simple PEC-modcl serves to not obscure the main results with too massive algebra. 
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3.1 Surface currents on a wide strip 



Here, we consider a comparatively wide strip with the half-width h = 0.5 m. The strip 
is at the height a = 1 m and the dipole source is at the height Zq = 0.5 m, above the 
PEC plane. The required number of spectral coefficients scales with h/X and here we use 
nT-max = 3 in (15^-(1571). i.e. the sought coefficients are {c„,(i„}^^Q. 

For graphical clarity and space reasons, we only consider the portion -3 m < x < 3 m of 
the strip, and the results shown are only for the real part (one particular instantaneous 
value) of the surface current density. 

In cases with the dipole moment in the x-direction, the field will be purely TM in the 
x-direction, whereby Ky = 0. This extra numerical check resulted in i^j^-values that were 
about five orders in magnitude smaller than the i^a,-values. Hence, for x-directed dipoles 
only results for are presented. 

We check the symmetries expected when the dipole is centered with the strip at yo = 0. 
The results are shown in Figure [21 We see that Kr,; and Ky exhibit the symmetries 



expected, due to the locations and orientations of the dipole sources. In Figure 2(b) 



we see that the y-directed dipole yields odd symmetry in K^, resulting in a zero net 
current fiowing along the strip (cf. f l59|) ). Also, Ky dominates over K^, since at the strip 
the exciting dipole electric field is predominantly in the ^/-direction. For the z-directed 
dipole, we see in Figure 2(c) that, except for in the region above the dipole, K^ dominates 
over Ky. 



3.2 Surface currents on a narrow strip 



Here, we consider a comparatively narrow strip with the half-width h = 0.1 m. The strip 
is at the height a = 6 m and the dipole source is at the height zq = 5.5 m, above the PEC 



plane, and off-centered at z/q = 0.5 m. Here we use mmax = 1 in flMj) - fl57j) . i.e. the sought 
coefficients are {c„, dn}n=o- 

The conjecture is that for a narrow strip the longitudinal! current dominates over the 
transverse current, i.e. \Kx\ ^ \Ky\. Since an x-directed dipole yields Ky = 0, results 
for that case are omitted, due to space reasons. The results for y- and ^-directed dipoles. 



shown in E] Figures 3(a) and 3(b) respectively, show that for both the y-directed and the 



2r-directed dipole, K^ dominates over Ky, which supports the conjecture that on a narrow 
strip the surface current density is predominantly longitudinal. 



3.3 Simplified model for narrow strips using only the total long- 
itudinal currents 

Guided by the results in subsection 13. 2[ we will here make approximations that result in 
a considerably simplified model. 

First, the transverse surface current is neglected completely, i.e. we set Ky = whereby 
{dn = 0}^Q (cf. (140!) ). Next, assuming that the strip is not too close to the PEC plane 

■^Note that "longitudinal" here refers to the direction of the strip, the cc-direction, while in Section [2] 
"longitudinal" refers to the decomposition direction of the fields, the z-dircction. 
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Figure 2: Surface current when h = 0.5 m, a = 1 m. Dipole at = 0.5z m. 



and not too close to the dipole, we neglect proximity effects and assume that in ( 139|) the 
longitudinal surface current density can be described by the zeroth order coefficient only, 
i.e. {c„ = 0}^p From (l39l) and (l59l) . it thus follows that 



and, introducing the spectral current I{k. 



I{x) H(/i-||/|) 



Tih 



nhco^kx), that fjlTj) reduces to 
i^,(fct) = /(A;,)Jo(A;,/i) 



(60) 



(61) 



In order to avoid an over-determined system, the final approximation is that we omit 



usmg 



to enforce -E^ = on the strip. Hence, we only use m = in 



to enforce 



that the longitudinal electric field component = on the strip, whereby (150|) reduces 
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-3-2-1 1 2 3 

(a) Dipole moment p ~ ly Cm. 




(b) Dipole moment p = Iz Cm. 

Figure 3: Surface current when /i = 0.1 m, a = 6 m. Dipole at vq = 0.5y + 5.52 m. 
Note that, for graphical clarity, different scalings are used on the x- and ?/-axes. 



to 



Ex{K, ky) Jo{kyh) dky = 



(62) 



Now, under these approximations, the general equations fl55l) - fl57|) fall out, leaving us with 
f l5^ . which, as a scalar equation for the spectral current, becomes 



I{k, 



e - ki r 1 - e-j2'=-'^ 



er] 



e - kl 



-]kz\a~zo\ 



-]kz(a+zo) 



Joi^kyh) dky < 



k 

kx k' 



] cos{kyyo) : p = px 



e — e 

^ ^g-jfc.|a-^o| _ g-jfc.(a+^o)j j sinikyVo) ■.p = py 



k^ [sgn(a — Zq) e j'"^'" + e j'=^("+^o)j cos{kyyQ) : p = pz 

(63) 



To verify the approximation of using only the total longitudinal current, we define the rel- 
ative difference d{x) = |/app(x) — /gcn(a;)| / \Igcn{x)\, where /app(x) is the current obtained 
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Figure 4: Relative difference between currents obtained using the approximate method 
and the general method. 

from using the approximate equation and Igen{x) is the current obtained from using 
the complete system of equations. In the comparison example, the geometrical parame- 
ters are the same as in Figure [31 With the dipole at x = 0, symmetry yields that we only 
need to consider a; > 0, and here we consider the range < a; < 40 m. The results are 
shown in Figure HI where we see that d{x) <C 1, even in this case with a strip as wide as 
1/5 of the wavelength. 

3.4 Results for narrow strips compared with results for thin 
circular wires 

In this section, we compare results for the narrow strip with results for the equivalent cir- 
cular wire. With negligible proximity effects, the per length inductances and capacitances 
for a wire with radius s become the same as for the strip with halfwidth h ii s = h/2; see 



3.4.1 Identification of TEM-mode in tiie x-direction 

For a wire/strip at small height (in terms of wavelengths) above the ground plane, the 
TEM-mode (if excited) is expected to be the dominating mode, in the field and cur- 
rent distributions. Note that the propagation along the wire/strip is in the x-direction, 
wherefore TEM here refers to field components in the y- and z-directions only. 

In Appendix [HI we derive the expressions f ll20l) - fll22p for the ideal TEM-mode currents 
on the wire excited by dipole sources in the principal directions. 

In Figures [MSI we present results for the excited currents with the strip at the heights 
a = 1.0 m and a = 2.0 m, respectively (the wavelength is 1 m). In the same figures, 
we also present the TEM-mode current on the equivalent circular wire, as well as the 
difference between the total current on the strip and the TEM-mode current on the wire. 

(11201) yields that x-directed dipoles will not excite TEM-modes on the wire. Correspond- 
ingly, we obtain in the expression f[63l) for the spectral current on the strip that for an 



e.g. m- 
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x-directed dipole there is no pole aX = k (the TEM-mode wavenumber) , whereby the 
TEM-mode is not excited. Hence, the results in Figures 5(a) and 6(a) are for higher order 
modes on the strip only. With a = 1.0 m we see in Figures 5(b) and 5(c) that the excited 
strip current is explained reasonably by the wire TEM-mode current, apart from in the 
vicinity of the dipole. Roughly, the difference between the total current on the strip and 
the wire TEM-mode current is the contribution from the higher order modes on the strip. 

Increasing further to a = 2.0 m we see in Figures 6(b) and 6(c) that, within the distances 
considered, the strip current is no longer dominated by the TEM-mode. 



3.4.2 Comparison with results obtained using MoM-solver 

Since the spectral method is designed to handle strips of infinite length, and commercial 
softwares that can handle infinitely long wires are not readily available, it is difficult to 
set up a proper comparison example. Nevertheless, we use the method of moment code 
NEC [31] to compare the results for a wire of finite length with the results using the 
spectral method. The finite length wire is located at |x| < 40 m and is left open-ended. 

In Figure [7] we see that with an x-directed dipole the agreement is quite satisfactory. The 
reason is the absence of the TEM-mode. Due to radiation, the higher order modes suffer 
from a much faster attenuation, resulting in quite small reflections from the open-ended 
wire in the NEC-model. On the other hand, in Figure |8] we see that with a 2;-directed 
dipole the agreement is far from satisfactory. Here, the TEM-mode is present, and due to 
its small attenuation the reflections from the open end-points of the wire creates standing 
waves in the NEC-model. Along the infinitely long strip there are no standing waves (but 
of course interference between the strip and the PEC plane). 



4 Discussion and conclusions 

In this paper, we have investigated high frequency wave propagation for the generic case 
of a metallic PEC strip above a PEC ground plane. Starting from a general formulation 
for a wide strip, we have verified the approximation of using only longitudinal currents 
on a narrow strip. Furthermore, we have verified the equivalence between narrow strips 
and thin wires, when proximity effects can be neglected. 

Although generic, the case studied includes the key-parts regarding the mathematical 
analysis and the approximations made. Hence, it is in principle straight-forward to gen- 
eralize to multiple parallel strips, multiple layers (including lossy media) and multiple 
dipole or distributed sources. However, when considering wide strips, requiring the com- 
plete system f lM|) -f lFr|l of equations, a generalization can be rather cumbersome. On the 
other hand, the simplified model for narrow strips (equivalent with thin wires), based 
on equation (!63|l . is more tractable under a generalization to multiple wires and/or lay- 
ers. Hence, the method outlined in this paper has potential usefulness for estimating 
long range propagation of high frequency waves in wire structures. Such structures can 
be power lines and railway feeding systems, containing both over-head wires and wires 
submerged into ground. 
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A Derivation of transversal field expressions 



Here, we derive the expressions Qj and ffTOj) for the transversal fields. 
With 

i;±(fct,^),/f^(fct,^)oce^j^^^ (64) 

and @ to express the fields in the source- free (J=0) Maxwell's equations ([1]) and ([2]), 
the spectral fields fulfil 

(fct ± k,z) X = kr]H^ (65) 

(fct ± k,z) X r]H^ = -kE^ (66) 

Decomposing fl65|l . f l66|l and the fields into longitudinal and transversal parts, we obtain 

kt X Ef = kj]Hfz (67) 

kt X zEf ± k^z X i;^^ = kriH^ (68) 

fct X r]Hf = -kEfz (69) 

fct X zr]H^ ± x r]Hf = -kEf (70) 



Performing the operations i:x f[68|) and zx fjTOj) . we obtain 



ktEf T k,E, = kzx r]Ht (71) 
ktrjH^ T k^riH, = -kz x E, (72) 



Eliminating z x E'^ from ( 1^ and (172]) . and using ([H]), we obtain OH]): 

Ef = ^ [±A;,fct^± - A;i; X k.r^Hf] , ©' 

Similarly, eliminating z x ?7i?^ from ( !70|) and (!7T|) . and using ([8]), we obtain ( 1101) : 

r/J/t^ = ^ [±k,k,r^Hf + kzx k,Ef] , ([lOD' 

Note that the outlined derivation is the common way of expressing the transversal fields 
in e.g. a waveguiding structure, where one instead has ki — )■ jVt (see e.g. [26], Section 
5.1). 
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B The spectral fields from the electric dipole 



In this appendix, we derive the expressions fl23|) - fl25|) for the functions Ai{kt,z) and 
Bi{kt, z), in the spectral fields from the electric dipole. 

B.l Spectral fields from a localised source 

First, we consider the more general case with a distributed but localized source J{r) 
located in the slab region between the planar surfaces Si, at z = zi, and S2, si z = Z2 
{zi < Z2)', see Figure ini From the losses in the medium and conservation of energy, it 
follows that the modes attenuate away from the slab region containing the source. 

Hence, from fll5p -( l20p it holds for z > Z2 that 



^ ^ [A+{k,) /+(fct) + B+{k,) ^(fct) ]e-j^^^e-j'=-''d'A:t (73) 

r]Hir) = -i?+(fct)/+(fet)]e-j'=^^e-j*'-''d2A;t (74) 

while for z < zi 

= A f [^"(fct)r(fct) + 5"(fct)^(fct)]e+j'=^^e-j'^-''d2A;t (75) 

rjHir) = J J- A~{k,) g{k,) + {k,) f- {k,)]e'^^'^%-^'^-f>d'h (76) 

From ([3]) and (jlj), we have the completeness relation: 

_L j g±j(fct-fcO Pds = 5(fct - K) (77) 

where 6() denotes the Dirac delta distribution. From (1771) . we have the corollaries 

[f^ih) X gik[)] ■ ^ ^C^-^O-Pd. = -M_k, . k[8ik, - k[) = ^Sih - K) 

(78) 

[/±(fcO X f^{k[)] ■ z-i^ f eiC^-^^^-^d. = (79) 

[/±(fcO X r (fcO] ■ I eiC^-^^^-^ds = (80) 

J s 

[g{K) X g{k[)] ■ I ei^^-^^^-^d. = (81) 

J s 

The Lorentz reciprocity theorem, see e.g. [32], yields 
j [Ex if*""* - E'''^' xH]-zds- f [Ex iJ*"'* - E""'' x H] ■ zds = f J ■ E'^'^'dv 



(82) 
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where ^E^'^^^, JJ^^^^^ is a solution to the Maxwell's equations that is source-free in the slab 
region between Si and S2- It has been assumed that the contribution to the total surface 
integral from the closing cylindrical boundary at infinity vanishes, for which a sufficient 
condition is that the source is localized and that the medium is lossy. 

Next, the mode coefficients will be calculated by inserting (173|) -( [76|) into (182|) and using 
the results (175]) -(IgT]). 

With = /+(fct)e-j'=^V*'* '',r/Jf*''* = g{kt) e-^''^'e^''' P, we obtain 



A+{k[) r{k[) X g{k,) + B+{k[) g{k[) x g{k,) 



ze 



A-{k[) f-{k[) X g{k,) + B-{k[) g{k[) x g{k,) 



JC 



A-{k[) /+(fct) X gik[) - B-{k[) rih) X f-ik[) 
-J^A-{k,) = 7^ J J{r) ■ f+{k,)e"^'^%^''-f>dv 



47r2 



^3) 



With E^^'^ = f-{kt)e+^''^'e^''' P,r]H^^''^ = -c/(fet) e+j^^V*'*''', we obtain 



/ d^A:; [ - A^{k[) r{k[) X g{k,) - B+{k[) g{k[) x g{k,) 
- A-^{k[) r (fcO X g{k[) + B^{k[) r{kt) X r{k[) 



d'k' 



K 



A-{k[) r{k[) X g{K) - B~{k[) g{k[) x g{k,) 



+ A-{k[) /-(fct) X g{k[) - B-{k[) r{kt) X f-{k[) 
= —g^A+{k,) = rij J(r) ■ /^(fct)e+j^^-V^-''dt; 



ze 



47r2 



With E''''' = g{kt)e-^^^'(^'''-P,riH 



test 



-/■^(fct)e j'^^^e^'^fP, we obtain 



d'k[ - A^k[) r (fco X r (fco - B^k[) g{k[) X r (fco 



- A+{k[) g{k,) X g{k[) + B^k[) g{k,) x r{k[) 



ze 



47r2 



dX -A-ik[)r{K) X rik,)-B-ik[)g{k[) X 



JC 



+ A-{k[)g{k,) X gik[) - B-{k[)g{k,) x r{K) 



gj(fct-fcO-Pds 



'>k k r 

'-^B-{k,)=r] J J(r)-^(fcOe-j'=^V^-''d^ 



^5) 
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With i;*"'* = g{kt)e+^''-'(^''^ f,7]H^'''^ = /"(fct) e+j'^-^eJ'^f'', we obtain 



- A^{k[) g{h) X g{k[) + B^{k[) g{k,) x 



4vr2 



d^A;; A-{k[)r{K) X r{k,) + B-{k[)g{k[) x r(fcO 



/c 



+ X - B~{k[)g{h) X r (fco 

y^B+{kt) = 7] J J(r) ■ ^(fct) e+j'^^ V^-^dw 



47^2 



Summarizing, and using (fT2l) and ( fT3i) . we obtain the relations 

A^(fct) = {k.k, ± . I J(r) ei^^-^^^^-^Mt 



^6) 



^7) 



5±(fct) = (z xkt)- j J(r) ej^'^'-P^'^^^My 



An electric dipole p at the location Vq = Pq + zqZ gives rise to the current density 

J{r)=]u:pb{r-r^) (89) 
where 6() denotes the Dirac delta distribution. 

Using f l89|) in f l87|) and flHHj) . we obtain for the three principal orientations of the dipole 



p = px 



p = py 




ipkx 

1 

2e 
2ek, 



(90) 



ipki 



2e 
YpkK 
' 2ek, * 



(91) 



p = pz 



JPkj UkfPr,±k,Zo) 



'2eK 



(92) 



Finally, with 

A(fct,^) = A±(fcOe^j^^^ 
we obtain the expressions (12^ - fl2S]) . 



B,{k,,z)=B^{k,)e 



^TikzZ 



(93) 
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C Numerical evaluation of spectral integrals 



Here, we describe the numerical evaluation of the spectral integrals in flMl) -fl5< 

Since the integrals run over the semi-infinite interval < ky < oo and are not open for 
easy analytic evaluation, they must be approximated over a finite interval < ky < ky^^^x- 
To increase the rate of convergence with increasing ky^^^x, we identify and subtract off 
leading behaviours of the integrands for large ky that can be integrated analytically. 

First, the longitudinal wavenumber is 



kz 



T^'2i Z-^ 2 



(94) 



With Im{A;^} < and Q < ky < oo, it holds for ky ^ \kx\ , \k\ that kz ~ —']ky. 

Hence, for large ky we have in the left-hand sides of (15^ - (!571) integrands that behave as 



1 - e 



-j2fcza 



ky. 



- Jj/j {kyK) Jjyj {kyh) 



1 — e"^''!'" J 

- J^i (kyh) Jj,2 (kyh) ^ — J^^ (kyh) Ji,2 (kyh) (95) 



k,. 



(l - e'j^*-'"") - -j 

— ^y^{kyK) 3u2{^y^) ~ ~j "^uiikyh^ 3i^2(^kyh^ 

Ky Ky 



(96) 



From [21], formula 6.574, we have the following results for non-negative integers m and 



n: 







•i2m{kyh) ^2n{kyh) 



[kyh) {kyh) 



dk 



dk 



^ 2{m + n)' 



(m,n)^(0,0) 



y 



.y 2 (m + n + 1) 

Using (jnZD, the integrals in (IMl) and (1551) are approximated as 



(97) 
(98) 







kz 



-hm{kyh) hn{kyh) dky 



I _ Q-i2k^a j 

/c? k^i 



J2m 

(/Cy/l) J2n{kyh) dky + 



j 6^ 



2 (m -|- n) 



(m,n)7^(0,0) (99) 



l^yl^Z 



(1 



hmikyh) J2n{kyh) dky 
{k' - k'y) , j 



ky 



hm{kyh) J2n{kyh) dky 



j 5r 



2 (m + n) 



(100) 



For the case m = n = 0, the first integral in (1541) is approximated as 



1 - e 



-i2kza 



l/h ^ _ g-j2fc,a 



k. 



1 - e 



oo t2 



Jo(«) 



dw, 



(101) 
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where / ^^dn ^ 0.3438831082, evaluated as a generalized hypergeometric function. 
Ji u 

Using ( p8i) . the integrals in (l56i) and (1571) are approximated as 



>2m+l 



(kyh) hn+i{kyh) dky 



-j2fcza 



yj 



'2m+l 



(kyh) J2n+i{kyh) dky + 



j 5r 



2 (m + n + 1 



(102) 



(1 _ e-j2fc.a^ (^2 _ ^2) 



^2m+i{kyh) J2n+i{kyh) dky 



(1 - e-j^^-") (fc^ - fcg) ^ j 



J2m+l(^y^) J2n+l(^y/i) d/Cy — 



j 5r 



2 (m + n + 1) 

(103) 



In the right hand sides of fl5^ - fl56l) . the exponential functions typically guaranty rapid 
convergence of the integrals, except for small values of a — zq, i.e. when the dipole source 
is close to the plane containing the strip. To account for such cases, we for x-directed 
dipoles extract the following explicit integral for the leading behaviours ( |2l] 6.621): 



•^mikyh) dky 



1 

m 



^ , m m + 1 




m 



1,2, 



(104) 



where F() is the hypergeometric function. For y- and ^-directed dipoles we extract the 
following explicit integral for the leading behaviours ( [21] 6.611): 



-k^b 



" '^mikyli) dky 



m 



0,1, 



(105) 



In both f ll04p and fll05p . the complex variable h has Re {6} = 6 =F and Im {6} = ±yo- 

The integrals in ( I54l) -( l57|) have been computed numerically by the routine "quadl" in the 
Matlab software. We notice that the integrands have poles where /c^ = 0. By considering 
the medium as generally lossy, the constitutive parameters have non-zero imaginary parts, 
moving the poles off the real axis. However, in the near lossless limit these poles are close 
to real axis, causing numerical problems in the integrals. To alleviate this problem, 
the numerical approximations of the integrals in (ISSj) are calculated by deforming the 
integration path to the contour shown in Figure [TUl 

With Ak^ ^ — Im{A;}, we obtain better numerical stability for the integrals in flMjl -f lFTj) . 
and for the subsequent integrals in (158|) . The upper limit on Ak^ is in practice determined 
by the acceptable exponential growth in the term ^^'^'^ (see (158|1 ). 
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D The TEM-mode current along an x-directed cir- 
cular wire 



In this appendix, we derive the expressions for the TEM-mode current on a circular wire, 
used in Section 13.4.11 for comparison purposes. 

Consider a PEC circular wire oriented in the x-direction. The wire has radius s and 
is centered at {y = 0,z = a) above the PEC ground plane (s < a). Note that in the 
present context, y and z are the transversal directions, subscripted with t, while x is the 
longitudinal direction. 

Applying f l82|) in the present context, the TEM-mode excitation coefficients become 



I'- 




Jv 




2 / {E, 


X Ht) ■ xda 


J A 



Dt = '-j^ (106) 



where V is the volume containing the source and A is the halfplane (z > 0)n(— oo < y < oo) 
expect for the cross-section of the wire. Et{y, z) and Ht{y, z) are the modal field patterns 
in the cross-sectional plane. 

From the TEM-mode relation = rj^^x x E^ (where a; ■ = 0), we obtain 



/ J-^;te±j^"dv 
Jv 

E^da 



= -l^r (107) 

J A 

Introducing, the scalar potential function $(?/,2), Et = — Vt$ and Vt ■ Et = V^^ = 
together with the Gauss' theorem in the plane yields 

/ E^da= I |Vt$|^da= / [Vt ■ ($Vt$) - <I>Vt $] da = * <I)n ■ Vt$d/ (108) 

J A J A J A Jc 

where C = Cq + Coo + Cg, in which Cq is the line z = that is closed with the curve Coo 
at infinity (with z > 0), and Cg is the circumference of the wire. With $ = on Co and 
Coo, with the constant value $ = $s on Cg, and the boundary condition a = enVt^ for 
the surface charge density, we obtain 

E^da = / -dl = (109) 

where A is the line charge density on the wire. 
By standard image theory |3l], we obtain 




Hy,z)^^Jn,l«4±4±4^4±^^^l^^ (110) 



27re y y2 _^ + a2 _ _ 2z\f - s 



On the wire, i.e. when y"^ + [z — a) = s^, flllOp yields 



^ A a + y/a^ — A /a\ 

$s = In = arccosh — 1 (111) 

zne s ZTTE \s, 
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From flllOp . we obtain 



E,{y, z) = -VMy, ^) = 7T^e(y, z) (112) 

zvre 

where 

_ yy+{z- a/q^ - g^) z _ + (2: + y/g^ - s^) 2: ^^^^^ 
?/2 + (2 — a? — s'^Y y'^ + (yZ + \/ a? — s^)^ 

With the dipole p at ro, yielding J(r) = ]ujpb{r — ro), it follows from f ll07p . f llOQp . fillip 
and (fml) that 



2Aarccosh(a/s) 
Hence, the //-field becomes 

To determine the current distribution along the wire, we must find the modal 
associated with 



(114) 



(115) 
current 



-X X Et{y,z) 
V 



A 



2neri 



yz — (^z — \J — s^) y yz — (^z + \/ — s^) ^ 
?/2 + (2; — a/o^ — s^)^ ?/2 _|_ (2 + \/ a? — s'^Y 



On the ground plane, the surface current density becomes 



K{y) = z^H,{y,z = Q) 



Tier] y"^ + a? — 



X 



which yields that the oppositely directed current on the wire becomes 

2A 



K{y) ■ xdy 



Trerj Jq y"^ + — 



dy 



A 



ST] 



(116) 



:ii7) 



:ii8) 



Hence, the TEM-mode current distribution on the wire becomes 



/tem(x) = -sgn(x - xo) 2Aarccosh(a/.) ^ 



2arccosh(a/ s) 



(119) 
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Explicitly, for the three principal directions of the dipole, it follows from (I113p and (11191) 
that 



p = px^ Item{x) = 



(120) 



p = py ^ 



Item{x) 



jujp sgn[x - xo) yo 



(121) 



2arccosh(a/ s) [^2 ^ (^^o - Va"^ - s^) ^ 



p = pz ^ 



jup sgn{x - Xo) zo - Va^ - s2 



zq + \J — 



(122) 



2arccosh(a/s) [^2 ^ _ ^^2 _ ^2^2 



vl + (2:0 + - s2) ^ 



The expressions (I12ip and (I122p hold for any wire radius < s < a, but if s ^ a, |a — 2;o|, 
we can use — ~ a and arccosh(a/s) ~ ln(2a/s), cf. (II lip . 
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Figure 5: Strip current wlien h = 0.02 m, a = 1 m. Dipole at ro = 0.5^ + 0.5z m. 
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(c) Dipole moment p = Iz Cm. 
Figure 6: Strip current when h = 0.02 m, a = 2 m. Dipole a.t vq = y + z m. 
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Figure 7: Strip current when h = 0.02 m, a = 6 m. Dipole at ro = Qy + 5.5z m. Dipole 
moment p = Ix Cm. Tlie reference results from tlie NEC-code are for an open-ended 
wire at Ixl < 40 m. 
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Figure 8: Strip current when h = 0.02 m, a = 6 m. Dipole at = 6y + 5.5z m. Dipole 
moment p = Iz Cm. Ttie reference results from tlie NEC-code are for an open-ended 
wire at \x\ < 40 m. 
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Figure 9: Geometry utilized when applying the Lorentz reciprocity theorem 
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Figure 10: Integration contour in the numerical approximation of 
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